#######################################################################
##### This file is part of the replication files of ###################
##### the article "Trading Integrity for Competence?" #################
##### , which draws Figure 1 and Online Appendix Figure A #############
#######################################################################

rm(list=ls())
library(readstata13)

library(ggplot2)
library(haven)
library(dplyr)
library(tidyverse)
library(magrittr)
library(lmtest)
library(sandwich)
library(pBrackets)
library(stringr)
library(gridExtra)
library(ggstance)
library(latex2exp)
devtools::install_github("eliocamp/tagger", force=TRUE)
library(tagger)

setwd("D:/Dropbox/Tradeoff replication")

data<-read.dta13("tradeoff_data.dta")%>%
  mutate(CI=(QK=="CIT"|QK=="CIM"),
         CC=(QK=="CCT"|QK=="CCM"),
         HI=(QK=="HIT"|QK=="HIM"),
         HC=(QK=="HCT"|QK=="HCM"))%>%
  mutate(municipal=rank==1,
         town=rank==0)%>%
  mutate(age=as.factor(age),
         income=as.factor(income),
         contact=as.factor(contact))


# As township is set as the reference, 
# CI, CC and HC coefs are township estimates (HI omitted)
est_town = "CI + CC + HC + municipal + CI:municipal + CC:municipal + HC:municipal"
est_ctrl_town = "CI + CC + HC + municipal + CI:municipal + CC:municipal + HC:municipal + male + han + age + hedu + income + ccp + urban + localcity +  hukouu  + gov + town_clean + town_corrupt + municipal_clean + municipal_corrupt  + contact"

## change reference group to municipality,
## so that CI, CC and HC coefs are municipal estimates (HI omitted)
est_municipal = "CI + CC + HC + town + CI:town + CC:town + HC:town"
est_ctrl_municipal = "CI + CC + HC + town + CI:town + CC:town + HC:town + male + han + age + hedu + income + ccp + urban + localcity +  hukouu  + gov + town_clean + town_corrupt + municipal_clean + municipal_corrupt + contact"


## Obtain estimates for Figure 1
output<-matrix(NA, nrow=9, ncol=6)%>%
  as_data_frame()  

output[, 1]<-rep(c(0, 1, 2), each=3)
output[, 2]<-rep(c("demotion", "promotion", "report"), 3)
output[, 3]<-0
output[, 4]<-0
output[, 5]<-0
output[, 6]<-c("(a)", "(b)", "(c)", "(d)", "(e)", "(f)", "(g)","(h)", "(i)")

for (i in 1:9){
  if (i<=3){ # township: (a) to (c)
    reg.town <-NULL
    reg.town <- lm(formula=paste(output[i,2], est_town, sep = "~"), data=data)
    output[i, 3]<- coeftest(reg.town, vcov = vcovHC(reg.town, type="HC1"))["CCTRUE", 1]
    output[i, 4]<- coeftest(reg.town, vcov = vcovHC(reg.town, type="HC1"))["CCTRUE", 2]    
  }else{
    if (i<=6) { # municipal: (d) to (f)
      reg.municipal <-NULL
      reg.municipal <- lm(formula=paste(output[i,2], est_municipal, sep = "~"), data=data)
      output[i, 3]<- coeftest(reg.municipal, vcov = vcovHC(reg.municipal, type="HC1"))["CCTRUE", 1]
      output[i, 4]<- coeftest(reg.municipal, vcov = vcovHC(reg.municipal, type="HC1"))["CCTRUE", 2]   
    }else{ # difference across levels: (g) to (i)
      reg.diff<-NULL
      reg.diff <- lm(formula=paste(output[i,2], est_town, sep = "~"), data=data)
      output[i, 3]<- coeftest(reg.diff, vcov = vcovHC(reg.diff, type="HC1"))["CCTRUE:municipalTRUE", 1]
      output[i, 4]<- coeftest(reg.diff, vcov = vcovHC(reg.diff, type="HC1"))["CCTRUE:municipalTRUE", 2] 
    }
  }
}

output<-output%>%
  set_colnames(c("rank", "outcome", "estimate", "se", "control",  "subfigure"))


## Include controls to get estimates for Online Appendix Figure A
output.ctrl<-matrix(NA, nrow=9, ncol=6)%>%
  as_data_frame()

output.ctrl[, 1]<-rep(c(0, 1, 2), each=3)
output.ctrl[, 2]<-rep(c("demotion", "promotion", "report"), 3)
output.ctrl[, 3]<-0
output.ctrl[, 4]<-0
output.ctrl[, 5]<-1
output.ctrl[, 6]<-c("(a)", "(b)", "(c)", "(d)", "(e)", "(f)", "(g)","(h)", "(i)")

for (i in 1:9){
  if (i<=3){ # township
    reg.town <-NULL
    reg.town <- lm(formula=paste(output.ctrl[i,2], est_ctrl_town, sep = "~"), data=data)
    output.ctrl[i, 3]<- coeftest(reg.town, vcov = vcovHC(reg.town, type="HC1"))["CCTRUE", 1]
    output.ctrl[i, 4]<- coeftest(reg.town, vcov = vcovHC(reg.town, type="HC1"))["CCTRUE", 2]    
  }else{
    if (i<=6) {
      reg.municipal <-NULL
      reg.municipal <- lm(formula=paste(output.ctrl[i,2], est_ctrl_municipal, sep = "~"), data=data)
      output.ctrl[i, 3]<- coeftest(reg.municipal, vcov = vcovHC(reg.municipal, type="HC1"))["CCTRUE", 1]
      output.ctrl[i, 4]<- coeftest(reg.municipal, vcov = vcovHC(reg.municipal, type="HC1"))["CCTRUE", 2]   
    }else{
      reg.diff<-NULL
      reg.diff <- lm(formula=paste(output.ctrl[i,2], est_ctrl_town, sep = "~"), data=data)
      output.ctrl[i, 3]<- coeftest(reg.diff, vcov = vcovHC(reg.diff, type="HC1"))["CCTRUE:municipalTRUE", 1]
      output.ctrl[i, 4]<- coeftest(reg.diff, vcov = vcovHC(reg.diff, type="HC1"))["CCTRUE:municipalTRUE", 2] 
    }
  }
}
output.ctrl<-output.ctrl%>%
  set_colnames(c("rank", "outcome", "estimate", "se", "control",  "subfigure"))


## Combine the estimates to draw Figure 1 and Online Appendix Figure A
output.all<-rbind(output, output.ctrl)

# Set new facet label names for rank (government levels)
rank.labs <- c("Township: \n CCT-HIT", "Municipal: \n CCM-HIM", "Difference: \n (CCM-HIM)-(CCT-HIT)")
names(rank.labs) <- c("0", "1", "2")

# Set new facet label names for outcome variable
outcome.labs <- c("Demotion", "Promotion", "Reporting")
names(outcome.labs) <- c("demotion", "promotion", "report")

## Figure 1: no controls
ggplot(data=output.all%>%filter(control==0), aes(x=estimate, y=1))+
  geom_point(size=2.5) +
  geom_linerange(aes(y= 1 , xmin=estimate-1.96*se, xmax=estimate+1.96*se), size=0.8) + 
  #geom_linerange(aes(y=1, xmin=estimate-1.64*se, xmax=estimate+1.64*se), size=1.2) + 
  geom_vline(xintercept = 0, lty=2, colour="grey", size=0.45) +
  geom_text(aes(x = estimate,  y = 1.25,  label = round(estimate, 3)), family="Times New Roman", size=3.5) +
  geom_text(aes(label = subfigure, family="Times New Roman", x= -Inf, y=Inf), hjust = 1, vjust = 1) +
  theme_bw() + xlab(NULL) +
  theme(legend.position = "bottom",
        panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        strip.text.x = element_text(size = 10),
        strip.text.y = element_text(size = 10),
        axis.text.x=element_text(size=10, colour="black", family="Times New Roman"), 
        axis.text.y=element_text(size=10, colour="black", family="Times New Roman",  angle=90, hjust=0.5), 
        axis.title.x=element_text(size=10, colour="black", family="Times New Roman", vjust=-1),
        axis.title.y=element_text(size=10, colour="black",family="Times New Roman"),
        text=element_text(size=10,  colour="black", family="Times New Roman")) + 
  scale_y_discrete(name ="", labels = "") +
  coord_cartesian(ylim = c(0, 2))+
  facet_grid(rank ~ outcome, scales = "free", 
             switch = "y", 
             labeller = labeller(rank = rank.labs, outcome = outcome.labs ))+
  tag_facets()

ggsave("Figure 1_without cotrls.jpeg", width = 15, height = 15, units = "cm", dpi=450)


## Online Appendix Figure A: with controls
ggplot(data=output.all%>%filter(control==1), aes(x=estimate, y=1))+
  geom_point(size=2.5) +
  geom_linerange(aes(y= 1 , xmin=estimate-1.96*se, xmax=estimate+1.96*se), size=0.8) + 
  geom_vline(xintercept = 0, lty=2, colour="grey", size=0.45) +
  geom_text(aes(x = estimate,  y = 1.35,  label = round(estimate, 3)), family="Times New Roman", size=3.5) +
  theme_bw() + xlab("") +
  theme(legend.position = "bottom",
        panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        strip.text.x = element_text(size = 10),
        strip.text.y = element_text(size = 10),
        axis.text.x=element_text(size=10, colour="black", family="Times New Roman"), 
        axis.text.y=element_text(size=10, colour="black", family="Times New Roman",  angle=90, hjust=0.5), 
        axis.title.x=element_text(size=10, colour="black", family="Times New Roman", vjust=-1),
        axis.title.y=element_text(size=10, colour="black",family="Times New Roman"),
        text=element_text(size=10,  colour="black", family="Times New Roman"))+
  scale_y_discrete(name ="", labels="")+
  coord_cartesian(ylim = c(0, 2))+
  facet_grid( rank ~ outcome, scales ="free", 
             switch ="y",
             labeller = labeller(rank = rank.labs, outcome = outcome.labs))+
  tag_facets()
ggsave("Figure A_with cotrls.jpeg", width = 15, height = 15, units = "cm", dpi=450)
